// standard headers
#include <iostream>
#include <ostream>
#include <istream>
#include <sstream>
#include <fstream>
#include <algorithm>
#include <map>
#include <vector>
#include <deque>
#include <string>

// boost headers
#include <boost/format.hpp>
#include <boost/tuple/tuple.hpp>

// component headers
#define _MY_STL_MACROS
#include <stlhelper.h>
#include <matrix.h>
#include <tree.h>
#include <tree_ops.h>
#define DEBUG_LEVEL 0
#include <debug_func.h>


using std::string;
using std::cerr;
using std::cout;
using std::vector;
using std::deque;
using std::ostream;
using std::ofstream;
using std::fill;
using std::istream;
using std::stringstream;
using boost::format;
using namespace bcpp_matrix;

namespace bcpp_tree
{
//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//	Operations

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888
struct do_map_UIDs_to_normalised_UIDs
{
	const bcpp_tree::t_tree&				tree;
	vector<unsigned>&						map_old_UIDs_to_new;
	unsigned								curr_UID;
	do_map_UIDs_to_normalised_UIDs(const bcpp_tree::t_tree& tree_, vector<unsigned>& map_old_UIDs_to_new_):
		tree(tree_), map_old_UIDs_to_new(map_old_UIDs_to_new_), curr_UID(0)
	{
		// allocate node number of UIDs
		map_old_UIDs_to_new.assign(tree.get_node_count(), 0);

		// assign UIDs for leaves first
		for (unsigned i = 0, sz = tree.get_leaf_count(); i < sz; ++i)
			map_old_UIDs_to_new[tree.get_leaf(i).get_UID()] = curr_UID++;
		recurse(tree.get_root());
	}
	void recurse(bcpp_tree::t_node* node);

};

void do_map_UIDs_to_normalised_UIDs::recurse(bcpp_tree::t_node* node)
{
	if (node->is_leaf())
		return;
	for (unsigned i = 0; i < node->get_children_count(); ++i)
		recurse(node->get_child(i));
	map_old_UIDs_to_new[node->get_UID()] = curr_UID++;
}

//________________________________________________________________________________________

//	do_sort_by_name

//________________________________________________________________________________________
struct do_sort_by_name
{
	void operator()(t_node* root) {recurse(root);}
	const std::string* recurse(t_node* node);
};

struct op_sort_by_str_ptr
{	bool operator()(const string* a, const string* b) const{return (*a < *b);} };


const string* do_sort_by_name::recurse(t_node* node)
{
	if (node->is_leaf())
		return &node->get_name();

	// only leaves don't have children
	assert(node->get_children_count());

	// get ultimate names of children (name = recursive depth first search to first leaf)
	typedef std::multimap<const string*, unsigned, op_sort_by_str_ptr> t_map_pstr_index;
	t_map_pstr_index names;
	for (unsigned i = 0; i < node->get_children_count(); ++i)
		names.insert(std::pair <const string*, unsigned>(recurse(node->get_child(i)), i));


	// sort children by ultimate name
	std::vector<t_node*> sorted_children;
	for (t_map_pstr_index::iterator i = names.begin(); i != names.end(); ++i)
		sorted_children.push_back(node->children[i->second]);

	assert(stlhelper::is_unique(node->children.begin(), node->children.end()));
	assert(stlhelper::is_unique(sorted_children.begin(), sorted_children.end()));
	sorted_children.swap(node->children);

	return names.begin()->first;

}






//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//	do_get_matrix

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//	do_get_matrix

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888
struct do_get_matrix
{
	do_get_matrix(t_dense_matrix_f & matrix_)
		:matrix(matrix_){}
	t_dense_matrix_f&				matrix;
	std::vector<float >	branch_lengths;
	void operator()(bcpp_tree::t_tree::t_node_type* root);
	void recurse(bcpp_tree::t_tree::t_node_type* node, std::deque<unsigned>& parent_branch);
};

void do_get_matrix::operator()(t_node* root)
{
	branch_lengths.assign(root->get_tree()->get_leaf_count(), 0.0);

	std::deque<unsigned> parent_branch;
	recurse(root, parent_branch);
}



void do_get_matrix::recurse(t_node* node, std::deque<unsigned>& parent_branch)
{
	assert(node);

	// save index in parent's list if leaf
	if (node->is_leaf())
	{
		assert((int)branch_lengths.size() > node->get_leaf_UID());
		parent_branch.push_back(node->get_leaf_UID());
		branch_lengths[node->get_leaf_UID()] += node->get_length();
		return;
	}

//	assert( node->get_children_count() == 2);

	// create lists of indices for all descendents of each child
	std::deque<std::deque<unsigned> >	branches(node->get_children_count());
	for (unsigned i = 0; i < node->get_children_count(); ++i)
		recurse(node->get_child(i), branches[i]);

	assert(branches.size());
	for (unsigned i = 0; i + 1 < branches.size(); ++i)
	{
		for (unsigned j = i + 1; j < branches.size(); ++j)
		{
			std::deque<unsigned>& indices1 = branches[i];
			std::deque<unsigned>& indices2 = branches[j];
			for (unsigned ii = 0; ii < indices1.size(); ++ii)
			{
				for (unsigned jj = 0; jj < indices2.size(); ++jj)
				{
					unsigned gene_1 = indices1[ii];
					unsigned gene_2 = indices2[jj];
					assert(gene_1 < matrix.size());
					assert(gene_2 < matrix.size());
					matrix[gene_1][gene_2] =
					matrix[gene_2][gene_1] =
						branch_lengths[gene_2] +
						branch_lengths[gene_1];
				}
			}
		}
	}

	// now move all into the branch of my parent
	for (unsigned i = 0; i < branches.size(); ++i)
		parent_branch.insert(parent_branch.begin(), branches[i].begin(), branches[i].end());
	for (unsigned i = 0; i < parent_branch.size(); ++i)
		branch_lengths[parent_branch[i]] += node->get_length();
	return;
}



//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//	t_node

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888
t_node::t_node(t_tree* tree_, unsigned UID_, int leaf_UID_, const string& name)
		:	tree(tree_), length(0.0),
			//len_to_root(0.0),
			//len_to_leaf(0.0),
			//level(0),
			//node_type(eNONE),
			parent(0), UID(UID_), leaf_UID(leaf_UID_)
{
	debug_func_cerr(5);
	set_name(name);
	debug_cerr(5, "UID=" << UID << ",Name=" << get_name());
}

void t_node::print(ostream& os) const
{
	os  << "NODE DUMP"
	<< "\n\tID    = "   << get_name()
//	<< "\n\ttype  = "   << enode_type_names[get_node_type()]
//	<< "\n\tlevel = "   << get_level()
	<< "\n\tlabel = "   << get_label()
//	<< "\n\tlength= "   << get_len_to_root()
	<< "\n";
}

void t_node::set_parent(t_node& new_parent)
{
	debug_func_cerr(5);
	// make sure parent is different
	if (&new_parent == parent)
		return;
#if (!defined(NDEBUG))
	debug_cerr(5, "UID=" << UID << ",Name=" << get_name());
	if (parent)
	{
		debug_cerr(5, "Reparent from " << parent->UID  << " to " << new_parent.UID);
	}
	else
	{
		debug_cerr(5, "Parent to " << new_parent.UID);
	}
#endif

	if (parent)
	{
		assert(includes(parent->children, this));
		parent->children.erase(find(beg_to_end(parent->children), this));
	}
	new_parent.children.push_back(this);
	parent = &new_parent;
}

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//	t_tree

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888
//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//	Operations

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888




//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888

//	tree definitions

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888



//________________________________________________________________________________________
//	c'tor
//________________________________________________________________________________________
t_tree::t_tree() :
					root(0) ,
					current_node(0),
					current_parent(0),
					finished(false),
					sort_by_names_alphabetically(true)
					//,max_level(0)
{
	add_node(false, "root");
	root			=
	current_node    =
	current_parent  = &shared_mem.back();
}


//________________________________________________________________________________________
//	join_nodes
//________________________________________________________________________________________
t_tree::t_node_type& t_tree::join_nodes(t_node_type& node1, t_node_type& node2)
{
	debug_func_cerr(5);
	t_node_type& new_parent = add_internal_node();
	reparent_to_grandparent();
	debug_cerr(5, "Merge: "
					<< node1.get_UID()
					<< " and "
					<< node2.get_UID()
					<< " with "
					<< new_parent.get_UID());


//	cerr <<
//			format("Merge: [%1$s] %|30t|[%2$s] => %|45t|[%3$s]\n")
//			%	node1.get_name()
//			%	node2.get_name()
//			%	new_parent.get_name();
	node1.set_parent(new_parent);
	node2.set_parent(new_parent);
	return new_parent;
}


//________________________________________________________________________________________
//	add_leaf
//________________________________________________________________________________________
t_tree::t_node_type& t_tree::add_node(bool is_leaf, const string& str, int leaf_UID)
{
	assert (!finished);
	assert (is_leaf || leaf_UID == -1);

	// do this BEFORE creating node as it might access these parameters
	//	in its constructor
	node_labels.push_back(string());
	node_names.push_back(string());


	if (is_leaf)
	{
		//t_node(tree_, UID, is_leaf, name)
		shared_mem.push_back(t_node_type(this, shared_mem.size(),
											// only set non-incremental
											// leaf_UID if specified
											(leaf_UID != -1) ?
											leaf_UID :
											leaves.size(),
											str));
		// save pointer to leaf
		leaves.push_back(&shared_mem.back());
	}
	else
		// leaf UID always == -1 if non-leaf node
		shared_mem.push_back(t_node_type(this, shared_mem.size(), -1, str));

	current_node    = &shared_mem.back();
	if (current_parent)
		current_node->set_parent(*current_parent);
	return *current_node;
}


t_tree::t_node_type& t_tree::add_leaf(const string& str, int leaf_UID)
{
	return add_node(true, str, leaf_UID);
}

t_tree::t_node_type&  t_tree::add_internal_node()
{
	// internal name
	stringstream os;
	os << "internal_node{" << shared_mem.size() << "}";

	add_node(false, os.str());
	current_parent = current_node;
	return *current_node;
}


//________________________________________________________________________________________
//	finish
//________________________________________________________________________________________
bool t_tree::finish()
{
	debug_func_cerr(5);
	finished = true;

	if (current_node != root && current_parent != 0)
		return false;

	if (sort_by_names_alphabetically)
	{
		do_sort_by_name sort_by_name;
		for_each(sort_by_name);
	}

	return true;
}

//________________________________________________________________________________________
//	label_current_node
//________________________________________________________________________________________
void t_tree::label_current_node(const string& str)
{
	debug_func_cerr(5);
	assert (!finished);
	if (current_node)
		current_node->set_label(str);
	else
	{
		if (current_parent)
			current_parent->set_label(str);
	}
}

//________________________________________________________________________________________
//	reparent_to_grandparent()
//________________________________________________________________________________________
bool t_tree::reparent_to_grandparent()
{
	debug_func_cerr(5);
	assert (!finished);

	// return false if already at root
	if (!current_parent)
		return false;
	assert(current_parent);
	current_node    = current_parent;
	current_parent  = current_parent->parent;
	return true;
}

//________________________________________________________________________________________
//	set_branch_len
//________________________________________________________________________________________
void t_tree::set_branch_len(const string& str)
{
	debug_func_cerr(5);
	assert (!finished);
	float len = atof(str.c_str());
	current_node->set_length(len);
	debug_cerr(5, "len=" << len);
}

//________________________________________________________________________________________
//	set_branch_len
//________________________________________________________________________________________
void t_tree::set_branch_len(float len)
{
	debug_func_cerr(5);
	assert (!finished);
	current_node->set_length(len);
	debug_cerr(5, "len=" << len);
}



//________________________________________________________________________________________
//	reset
//________________________________________________________________________________________
void t_tree::reset()
{
	debug_func_cerr(5);
	root			=
	current_node    =
	current_parent  = 0;
	node_names.clear();
	node_labels.clear();
	shared_mem.clear();
	leaves.clear();
	finished		= false;
	add_node(false, "root");
	root			=
	current_node    =
	current_parent  = &shared_mem.back();
}


//________________________________________________________________________________________
//	to_matrix
//________________________________________________________________________________________
void t_tree::to_matrix(bcpp_matrix::t_matrix_f& data) const
{
	unsigned size = get_leaf_count();
	data.matrix.allocate(size, data.null_value);
	data.names = get_leaf_names();
	assert(data.matrix.size() == data.names.size());
	do_get_matrix get_matrix(data.matrix);
	for_each(get_matrix);
	for (unsigned i = 0; i < size; ++i)
		data.matrix[i][i] = 0.0;
}



//________________________________________________________________________________________
//	get_leaf_names
//________________________________________________________________________________________
std::vector<std::string> t_tree::get_leaf_names () const
{
	std::vector<std::string> names;
	for (unsigned i =0; i < leaves.size(); ++i)
		names.push_back(leaves[i]->get_name());
	return names;
}



void t_tree::assert_consistency(t_node* node) const
{
	assert(node);
	if (node->get_children_count() == 0)
		return;
	for (unsigned i = 0; i < node->get_children_count(); ++i)
	{
		if (node->get_child(i)->parent != node)
			throw t_invalid_tree("The parent field of node " + node->get_child(i)->get_name() +
								" does not point to its proper parent " + node->get_name() + "\n");
		assert_consistency(node->get_child(i));
	}
}

void t_tree::assert_consistency() const
{
	assert_consistency(root);
}

bool t_tree::reroot(unsigned child_index, unsigned grandchild_index,
					float& child_length, float& grandchild_length, float& other_length)
{
	if (root->get_children_count() < 2 ||
		root->get_child(child_index)->get_children_count() < 2)
	{
		return false;
	}
	unsigned other_index = (child_index+1) % 2;
	t_node_type* child = root->children[child_index];
	t_node_type* grandchild = child->children[grandchild_index];
	t_node_type* other = root->children[other_index];

	// save lengths
	other_length = other->get_length();
	grandchild_length = grandchild->get_length();
	child_length = child->get_length();

	// set lengths
	other->set_length(child->get_length() + other->get_length());
	child->set_length(grandchild->get_length() * 0.5);
	grandchild->set_length(child->get_length());
	
	// swap parents
	std::swap(other->parent, grandchild->parent);
	// swap children
	std::swap(root->children[other_index], child->children[grandchild_index]);

#if (!defined(NDEBUG))
	assert_consistency();
#endif

	return true;
}

bool t_tree::undo_reroot(unsigned child_index, unsigned grandchild_index,
					float child_length, float grandchild_length, float other_length)
{
	if (child_index			> root->get_children_count() ||
		grandchild_index	> root->get_child(child_index)->get_children_count())
	{
		return false;
	}
	unsigned other_index = (child_index+1) % 2;
	t_node_type* child	 = root->children[child_index];

	// swap children
	std::swap(root->children[other_index], child->children[grandchild_index]);
	// swap parents
	std::swap(root->children[other_index]->parent, child->children[grandchild_index]->parent);

	// set lengths
	t_node_type* other = root->children[other_index];
	t_node_type* grandchild = child->children[grandchild_index];
	child	 = root->children[child_index];
	other->set_length(other_length);
	child->set_length(child_length);
	grandchild->set_length(grandchild_length);

#if (!defined(NDEBUG))
	assert_consistency();
#endif

	return true;
}

void t_tree::reroot(t_node* new_root_child1, t_node* new_root_child2)
{
	cerr
		<< "Reroot between "
		<< new_root_child1->get_name() << " and "
		<< new_root_child2->get_name() << "\n";
	if (new_root_child1 == root ||
		new_root_child2 == root ||
		(new_root_child1->parent==root && new_root_child2 ->parent== root))
		return;
	cerr << "...\n";
	//
	// make sure new_root_child1 and new_root_child2 are the right way around
	//	i.e new_root_child2 is child of new_root_child1
	//
	if (new_root_child1->get_parent() == new_root_child2)
	{
		assert(new_root_child2->get_parent() != new_root_child1);
		std::swap(new_root_child1, new_root_child2);
	}
	if (new_root_child2->get_parent() != new_root_child1)
	{
		cerr
			<< "new_root_child1=" << new_root_child1->get_name() << "\n"
			<< "new_root_child2=" << new_root_child2->get_name() << "\n"
			<< "new_root_child2 parent=" << new_root_child2->parent->get_name() << "\n";
	}
	assert(new_root_child2->get_parent() == new_root_child1);

	// set the children of the new root to both point to root
	new_root_child2->parent = root;

	// set the child rather than parent to root initially, will be swapped in loop next
	new_root_child1->change_child(new_root_child2, root);
	t_node* node = new_root_child1;

	t_node* child_to_change = root;
	while (node != root)
	{
		// save next node to change
		t_node* next_node = node->get_parent();

		// swap child and parent of node
		node->change_child(child_to_change, node->parent);
		node->parent = child_to_change;

		child_to_change = node;
		node = next_node;
	}

	// now child_to_change = the child of old root
	//			which is along the path to the new root
	// find other child of original root and link them up
	t_node* old_root_other_child = root->get_sibling_of(child_to_change);
	old_root_other_child->parent = child_to_change;
	child_to_change->change_child(root, old_root_other_child);

	// now set children of new root
	root->children[0] = new_root_child1;
	root->children[1] = new_root_child2;
}





t_node* node_from_old_UID(unsigned oldUID, t_tree& new_tree, const std::vector<unsigned>& old_UIDs_to_new)
{
	debug_func_cerr(5);
	assert(oldUID < old_UIDs_to_new.size());
	assert(new_tree.shared_mem.size() > old_UIDs_to_new[oldUID]);
	return &new_tree.shared_mem[old_UIDs_to_new[oldUID]];
}

void t_tree::recursively_copy_tree(const t_node& node, t_node& new_node, t_tree& new_tree,
											const std::vector<unsigned>& old_UIDs_to_new) const
{
	debug_func_cerr(5);
	// set parent
	if (node.parent)
		new_node.parent = node_from_old_UID(node.parent->UID, new_tree, old_UIDs_to_new);

	// set children
	assert(new_node.children.size() == node.children.size());
	for (unsigned i = 0; i < node.children.size(); ++i)
		new_node.children[i] = node_from_old_UID(node.children[i]->UID, new_tree, old_UIDs_to_new);

	// call recursively
	for (unsigned i = 0; i < node.children.size(); ++i)
		recursively_copy_tree(*node.children[i], *new_node.children[i], new_tree, old_UIDs_to_new);
}

void t_tree::get_normalized_tree(t_tree& new_tree) const
{
	debug_func_cerr(5);
	assert (finished);
	new_tree = *this;

	// get normalised UIDs map and indices for reordering
	std::vector<unsigned> old_UIDs_to_new;
	do_map_UIDs_to_normalised_UIDs prepare_UIDs(*this, old_UIDs_to_new);
	std::vector<unsigned> indices(old_UIDs_to_new.size());
	for (unsigned i = 0; i < indices.size(); ++i)
		indices[old_UIDs_to_new[i]] = i;

	vector<string> old_names (new_tree.shared_mem.size());
	for (unsigned i = 0; i < new_tree.shared_mem.size(); ++i)
		old_names[i] = new_tree.shared_mem[i].get_name();


	// order data to new UIDs
	for (unsigned i = 0; i < new_tree.shared_mem.size(); ++i)
		new_tree.shared_mem[i].tree = &new_tree;
	order_by_sorted_indices(old_names,  indices);
	order_by_sorted_indices(new_tree.shared_mem,  indices);
	order_by_sorted_indices(new_tree.node_labels, indices);
	order_by_sorted_indices(new_tree.node_names,  indices);

	// set UIDs correctly
	for (unsigned i = 0; i < get_node_count(); ++i)
		new_tree.shared_mem[i].UID = i;

	// leaves point to the right places
	unsigned cnt_leaves = new_tree.get_leaf_count();
	new_tree.leaves.clear();
	for (unsigned i = 0; i < cnt_leaves; ++i)
		// save pointer to leaf
		new_tree.leaves.push_back(&new_tree.shared_mem[i]);
	//leaf_UID should already be correct

	// set nodes
	new_tree.current_node	=
	new_tree.current_parent	=
	new_tree.root			= node_from_old_UID(root->get_UID(), new_tree, old_UIDs_to_new);
	recursively_copy_tree(*root, *new_tree.root, new_tree, old_UIDs_to_new);
	#ifndef NDEBUG
	assert_consistency();
	new_tree.assert_consistency();
	#endif
}

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888


//	to_nh_str


//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888
std::string t_tree::to_nh_str() const
{
	std::ostringstream ostrm;
	bcpp_tree_ops::do_newhampshire_format_str nh_str(ostrm);
	for_each(nh_str);
	return ostrm.str();

}

//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888


//	load_tree_from_newhampshire_format_str


//8888888888888888888888888888888888888888888888888888888888888888888888888888888888888888
boost::tuple<bool, string>
			load_tree_from_newhampshire_format_str(const string& tree_str, t_tree& tree)
{
	debug_func_cerr(5);
	string error_str;
	string::const_iterator phylip_iter = tree_str.begin();
	string::const_iterator phylip_end = tree_str.end();
	string str;
	char lastevent = '\0';
	bool first_parenthesis_read = false;
	bool last_parenthesis_read = false;
	while (phylip_iter != phylip_end)
	{
		switch (*phylip_iter)
		{
			case '(':
				debug_cerr(5, "Event=" << *phylip_iter);
				lastevent = *phylip_iter;
				if (!first_parenthesis_read)
				{
					debug_cerr(5, "Read first parenthesis" << *phylip_iter);
					first_parenthesis_read = true;
					break;
				}
				str="";
				tree.add_internal_node();
				break;
			case ')':
			case ',':
				debug_cerr(5, "Event=" << *phylip_iter);
				if ( str.length())
				{
					if ( lastevent == ':' )
					{
						debug_cerr(5, "Set branch length = " << str << "\n");
						tree.set_branch_len(str);
					}
					else
					{
						debug_cerr(5, "add leaf: [ " << str << " ]\n");
						tree.add_leaf(str);
					}
				}
				//else
					//tree.add_leaf("");

				if (*phylip_iter == ')')
				{
					if (!tree.reparent_to_grandparent())
					{
						if (last_parenthesis_read)
						{
							debug_cerr(5, "Too many ')'");
							error_str = "mismatching parentheses: too many closing brackets";
							return boost::make_tuple(false, error_str);
						}
						else
						{
							debug_cerr(5, "Swallow last matching parenthesis");
							last_parenthesis_read = true;
						}
					}
				}
				str = "";
				lastevent = *phylip_iter;
				break;
			case ':':
				debug_cerr(5, "Event=" << *phylip_iter);
				if (lastevent == '(' || lastevent == ',' )
					tree.add_leaf(str);
				else if (str.length())
					tree.label_current_node(str);
				str = "";
				lastevent = *phylip_iter;
				break;
			case '\n':
			case '\r':
			case '\t':
			case ' ':
				break;
			case ';':
				if (!tree.finish())
				{
					debug_cerr(5, "Not at root" << *phylip_iter);
					error_str = "missing closing parentheses";
					return boost::make_tuple(false, error_str);
				}
				if (lastevent != ')' || str.length())
				{
					if (!str.length())
					{
						str += lastevent;
					}
					debug_cerr(5, "Must end with closing parenthesis" << *phylip_iter);
					error_str = "letters (\"" + str + "\") after closing parentheses";
					return boost::make_tuple(false, error_str);
				}
				return boost::make_tuple(true, error_str);
			default:
				if (!first_parenthesis_read)
				{
					error_str = "missing first parenthesis";
					return boost::make_tuple(false, error_str);
				}
				str += *phylip_iter;
		}
		++phylip_iter;
	}
	if (tree.finish())
	{
		error_str = "the closing semi-colon is missing";
		return boost::make_tuple(true, error_str);
	}
	error_str = "the tree ended prematurely";
	return boost::make_tuple(false, error_str);
}

//________________________________________________________________________________________

//	read_tree_from_stream

//		Will not allow more than one tree on a line

//________________________________________________________________________________________
//bool read_tree_from_stream(istream& istrm, bcpp_tree::t_tree& tree, unsigned& cnt_lines)
//{
//	debug_func_cerr(5);
//	string tree_name;
//	while (istrm)
//	{
//		// ignore empty lines
//		if (istrm.peek() == '\n')
//		{
//			istrm.ignore(1);
//			++cnt_lines;
//			continue;
//		}
//
//		// if line begins with >, that is tree name.
//		// N.B. First N -1 lines with > will be ignore (overridden by last instance)
//		if (istrm.peek() == '>')
//		{
//			getline(istrm, tree_name);
//			++cnt_lines;
//			tree_name.erase(0, 1);
//			tree.set_name(tree_name);
//			continue;
//		}
//
//		string tree_str;
//		if (getline(istrm,  tree_str, ';'))
//		{
//	    	bool success;
//			string error_str;
//			boost::tie(success, error_str) =
//							load_tree_from_newhampshire_format_str(tree_str, tree);
//			if (success)
//			{
//				cnt_lines += std::count(tree_str.begin(),  tree_str.end(),  '\n');
//				return true;
//			}
//
//			throw t_tree_parsing_from_stream_error(tree_name, error_str, cnt_lines);
//		}
//
//	}
//	return false;
//}


//________________________________________________________________________________________

//	read_tree_from_stream

//		Will not allow more than one tree on a line

//________________________________________________________________________________________
bool read_tree_from_stream(istream& istrm, bcpp_tree::t_tree& tree, unsigned& cnt_lines)
{
	debug_func_cerr(5);
	string tree_str;
	string tree_name;
	if (!read_treestr_from_stream(istrm, tree_str, tree_name, cnt_lines))
		return false;

	tree.set_name(tree_name);

	bool success;
	string error_str;
	boost::tie(success, error_str) =
					load_tree_from_newhampshire_format_str(tree_str, tree);
	if (!success)
		throw t_tree_parsing_from_stream_error(tree_name, error_str, cnt_lines);

	cnt_lines += std::count(tree_str.begin(),  tree_str.end(),  '\n');
	return true;
}

bool read_treestr_from_stream(istream& istrm, string& tree_str, string& tree_name, unsigned& cnt_lines)
{
	debug_func_cerr(5);
	tree_str =
	tree_name = string();
	while (istrm)
	{
		// ignore empty lines
		if (istrm.peek() == '\n')
		{
			istrm.ignore(1);
			++cnt_lines;
			continue;
		}

		// if line begins with >, that is tree name.
		// N.B. First N -1 lines with > will be ignore (overridden by last instance)
		if (istrm.peek() == '>')
		{
			getline(istrm, tree_name);
			tree_name.erase(0, 1);
			++cnt_lines;
			continue;
		}

		if (getline(istrm,  tree_str, ';'))
		{
			// nothing read
			if (!tree_str.length())
				return false;
			return true;
		}

	}
	return false;
}


std::ostream& operator <<(std::ostream& ostrm, const bcpp_tree::t_tree& tree)
{
	if  (!tree.get_name().empty())
		ostrm << ">" << tree.get_name() << "\n";
	bcpp_tree_ops::do_newhampshire_format_str nh_str(ostrm);
	tree.for_each(nh_str);
	return ostrm;
}

std::ostream& operator <<(std::ostream& ostrm, const bcpp_tree::t_node& node)
{
	bcpp_tree_ops::do_newhampshire_format_str subtree(ostrm);
	subtree(&node);
	return ostrm;
}

}// namespace



